{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Cellular automata\n",
    "\n",
    "Code examples from [Think Complexity, 2nd edition](http://greenteapress.com/wp/complexity2), Chapter 5\n",
    "\n",
    "Copyright 2016 Allen Downey, [MIT License](http://opensource.org/licenses/MIT)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "from __future__ import print_function, division\n",
    "\n",
    "%matplotlib inline\n",
    "\n",
    "import numpy as np\n",
    "import matplotlib.pyplot as plt\n",
    "\n",
    "import thinkplot"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "from thinkstats2 import RandomSeed\n",
    "RandomSeed(17)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Zero-dimensional CA"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Here's a simple implementation of the 0-D CA I mentioned in the book, with one cell."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "n = 10\n",
    "x = np.zeros(n)\n",
    "print(x)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "To get the state of the cell in the next time step, we increment the current state mod 2."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "x[1] = (x[0] + 1) % 2\n",
    "x[1]"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Filling in the rest of the array."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "for i in range(2, 10):\n",
    "    x[i] = (x[i-1] + 1) % 2\n",
    "    \n",
    "print(x)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "So the behavior of this CA is simple: it blinks."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## One-dimensional CA"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Just as we used a 1-D array to show the state of a single cell over time, we'll use a 2-D array to show the state of a 1-D CA over time, with one column per cell and one row per timestep."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "rows = 5\n",
    "cols = 11\n",
    "array = np.zeros((rows, cols), dtype=np.int8)\n",
    "array[0, 5] = 1\n",
    "print(array)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "To plot the array I use `plt.imshow`"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "def plot_ca(array):\n",
    "    cmap = plt.get_cmap('Blues')\n",
    "    plt.imshow(array, interpolation='nearest', cmap=cmap)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Here's what it looks like after we initialize the first row."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "plot_ca(array)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "And here's the function that fills in the next row.  The rule for this CA is to take the sum of a cell and its two neighbors mod 2."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "def step(array, i):\n",
    "    rows, cols = array.shape\n",
    "    for j in range(1, cols):\n",
    "        array[i, j] = sum(array[i-1, j-1:j+2]) % 2"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Here's the second row."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "step(array, 1)\n",
    "plot_ca(array)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "And here's what it looks like with the rest of the cells filled in."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "for i in range(1, rows):\n",
    "    step(array, i)\n",
    "\n",
    "plot_ca(array)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "For a simple set of rules, the behavior is more interesting than you might expect."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**Exercise:** Modify this code to increase the number of rows and columns and see what this CA does after more time steps."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Cross correlation"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We can step the CA through time more quickly using \"cross correlation\".  To see how it works, the first step is to replace the slice operator with array multiplication.\n",
    "\n",
    "This window selects the first three elements of an array:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "window = np.zeros(cols, dtype=np.int8)\n",
    "window[:3] = 1\n",
    "print(window)\n",
    "print(array[4])\n",
    "print(window * array[4])"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Then we can use `sum` and the modulus operator to compute the state of the first cell during the next timestep."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "sum(window * array[4]) % 2"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "To compute the state of the next cell, we shift the window to the right."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "window = np.roll(window, 1)\n",
    "print(window)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "And repeat the multiply-sum-modulus operations."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "sum(window * array[4]) % 2"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Now we can rewrite `step` using these operations."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "def step2(array, i):\n",
    "    rows, cols = array.shape\n",
    "    window = np.zeros(cols)\n",
    "    window[:3] = 1\n",
    "    for j in range(1, cols):\n",
    "        array[i, j] = sum(window * array[i-1]) % 2\n",
    "        window = np.roll(window, 1)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "And we can confirm that we get the same result."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "for i in range(1, rows):\n",
    "    step2(array, i)\n",
    "\n",
    "plot_ca(array)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "That sequence of operations is called a \"sliding dot product\" or \"cross correlation\", and NumPy provides a function that computes it.  So we can replace the `for` loop with `np.correlate`.  The parameter `mode='same'` means that the result has the same length as `array[i]`. "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "def step3(array, i):\n",
    "    window = np.array([1, 1, 1])\n",
    "    array[i] = np.correlate(array[i-1], window, mode='same') % 2"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "And the result is the same."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "for i in range(1, rows):\n",
    "    step3(array, i)\n",
    "\n",
    "plot_ca(array)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "So that's good enough for a CA that only depends on the total number of \"on\" cells, but for more general CAs, we need a table that maps from the configuration of the neighborhood to the future state of the center cell.\n",
    "\n",
    "The following function makes the table by interpreting the Rule number in binary."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "def make_table(rule):\n",
    "    \"\"\"Make the table for a given CA rule.\n",
    "    \n",
    "    rule: int 0-255\n",
    "    \n",
    "    returns: array of 8 0s and 1s\n",
    "    \"\"\"\n",
    "    rule = np.array([rule], dtype=np.uint8)\n",
    "    table = np.unpackbits(rule)[::-1]\n",
    "    return table"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Here's what it looks like as an array:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "table = make_table(150)\n",
    "print(table)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "If we correlate the row with the window `[4, 2, 1]`, it treats each neighborhood as a binary number between 000 and 111."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "window = [4, 2, 1]\n",
    "corr = np.correlate(array[0], window, mode='same')\n",
    "print(array[0])\n",
    "print(corr)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Now we can use the result from `np.correlate` as an index into the table; the result is the next row of the array."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "array[1] = table[corr]\n",
    "print(array[1])"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We can wrap up that code in a function:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "def step4(array, i):\n",
    "    window = np.array([4, 2, 1])\n",
    "    corr = np.correlate(array[i-1], window, mode='same')\n",
    "    array[i] = table[corr]"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "And test it again."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "for i in range(1, rows):\n",
    "    step4(array, i)\n",
    "\n",
    "plot_ca(array)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "How did I know that Rule 150 is the same as the previous CA?  I wrote out the table and converted it to binary."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "collapsed": true
   },
   "source": [
    "## The Cell1D object"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "`Cell1D.py` provides a `Cell1D` class that encapsulates the code from the previous section.\n",
    "\n",
    "Here's an example that runs a Rule 50 CA for 10 steps."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "from Cell1D import Cell1D, Cell1DViewer"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "rule = 50\n",
    "n = 10\n",
    "ca = Cell1D(rule, n)\n",
    "ca.start_single()\n",
    "ca.loop(n-1)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We can display the results:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "viewer = Cell1DViewer(ca)\n",
    "viewer.draw()\n",
    "\n",
    "plt.savefig('chap05-1.pdf')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Here's the Rule 50 table."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "print(ca.table)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Another example:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "rule = 150\n",
    "n = 5\n",
    "ca = Cell1D(rule, n)\n",
    "ca.start_single()\n",
    "ca.loop(n-1)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "viewer = Cell1DViewer(ca)\n",
    "viewer.draw()\n",
    "\n",
    "plt.savefig('chap05-2.pdf')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "And one more example showing recursive structure."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "rule = 18\n",
    "n = 64\n",
    "ca = Cell1D(rule, n)\n",
    "ca.start_single()\n",
    "ca.loop(n-1)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "viewer = Cell1DViewer(ca)\n",
    "viewer.draw()\n",
    "\n",
    "plt.savefig('chap05-3.pdf')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Rule 30 generates a sequence of bits that is indistinguishable from random:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "rule = 30\n",
    "n = 100\n",
    "ca = Cell1D(rule, n)\n",
    "ca.start_single()\n",
    "ca.loop(n-1)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "viewer = Cell1DViewer(ca)\n",
    "viewer.draw()\n",
    "\n",
    "plt.savefig('chap05-4.pdf')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "And Rule 110 is Turing complete!"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "rule = 110\n",
    "n = 100\n",
    "ca = Cell1D(rule, n)\n",
    "ca.start_single()\n",
    "ca.loop(n-1)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "viewer = Cell1DViewer(ca)\n",
    "viewer.draw()\n",
    "\n",
    "plt.savefig('chap05-5.pdf')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Heres a longer run that has some spaceships."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "rule = 110\n",
    "n = 600\n",
    "ca = Cell1D(rule, n)\n",
    "ca.start_random()\n",
    "ca.loop(n-1)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "viewer = Cell1DViewer(ca)\n",
    "viewer.draw()\n",
    "\n",
    "plt.savefig('chap05-6.pdf')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Exercises"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**Exercise:** This exercise asks you to experiment with Rule 110 and see how\n",
    "many spaceships you can find.\n",
    "\n",
    "1. Read the [Wikipedia page about Rule 110](https://en.wikipedia.org/wiki/Rule_110), which describes its background pattern and spaceships.\n",
    "\n",
    "2. Create a Rule 110 CA with an initial condition that yields the\n",
    "  stable background pattern.  Note that the CA class provides\n",
    "`start_string`, which allow you to initialize the state of\n",
    "the array using a string of `1`s and `0`s.\n",
    "\n",
    "3. Modify the initial condition by adding different patterns in the\n",
    "  center of the row and see which ones yield spaceships.  You might\n",
    "  want to enumerate all possible patterns of $n$ bits, for some\n",
    "  reasonable value of $n$.  For each spaceship, can you find the\n",
    "  period and rate of translation?  What is the biggest spaceship you\n",
    "  can find?\n",
    "\n",
    "4. What happens when spaceships collide?"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Solution\n",
    "\n",
    "# The following function makes a CA with the given rule and runs it for `n` steps.\n",
    "\n",
    "def run_ca(init, n=None, rule=110):\n",
    "    m = len(init)\n",
    "    n = m if n is None else n\n",
    "    \n",
    "    ca = Cell1D(rule, n, m)\n",
    "    ca.start_string(init)\n",
    "    ca.loop(n-1)\n",
    "    return ca"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Solution\n",
    "\n",
    "# Here's the background pattern.  Notice that this implementation doesn't get the borders quite right.\n",
    "\n",
    "from Cell1D import Wrap1D\n",
    "\n",
    "background = '00010011011111'\n",
    "\n",
    "ca = run_ca(background * 5)\n",
    "viewer = Cell1DViewer(ca)\n",
    "viewer.draw()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Solution\n",
    "\n",
    "# Here's a spaceship that translates right.  The parameters to `viewer.draw` trim off the borders.\n",
    "\n",
    "ship1 = '0001110111'\n",
    "\n",
    "ca = run_ca(background + ship1 + background * 3)\n",
    "viewer = Cell1DViewer(ca)\n",
    "viewer.draw(start=4, end=-1)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Solution\n",
    "\n",
    "# Here's one that translates left:\n",
    "\n",
    "ship2 = '1001111'\n",
    "\n",
    "ca = run_ca(background * 5 + ship2 + background * 3)\n",
    "viewer = Cell1DViewer(ca)\n",
    "viewer.draw(start=4, end=-1)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Solution\n",
    "\n",
    "# And here's one that stands still.\n",
    "\n",
    "ship3 = '111'\n",
    "\n",
    "ca = run_ca(background * 2 + ship3 + background * 2)\n",
    "viewer = Cell1DViewer(ca)\n",
    "viewer.draw(start=4, end=-1)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Solution\n",
    "\n",
    "# When these two ships collide, they pass through each other:\n",
    "\n",
    "init = background*4 + ship1 + background*3 + ship2 + background*3\n",
    "\n",
    "ca = run_ca(init, n=180)\n",
    "viewer = Cell1DViewer(ca)\n",
    "viewer.draw(start=4, end=-1)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Solution\n",
    "\n",
    "# When `ship1` hits `ship3`, it creates one `ship2` and one other spaceship we have not seen.\n",
    "\n",
    "init = background*3 + ship1 + background + ship3 + background*3\n",
    "\n",
    "ca = run_ca(init, n=120)\n",
    "viewer = Cell1DViewer(ca)\n",
    "viewer.draw(start=4, end=-1)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**Exercise:** The goal of this exercise is to implement a Turing machine.\n",
    "\n",
    "1. Read about Turing machines at http://en.wikipedia.org/wiki/Turing_machine.\n",
    "\n",
    "2. Write a class called `Turing` that implements a Turing machine.  For the action table, use the rules for a 3-state busy beaver.\n",
    "\n",
    "3. Write a class named `TuringDrawer` that generates an image that represents the state of the tape and the position and state of the head.  For one example of what that might look like, see http://mathworld.wolfram.com/TuringMachine.html.\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Solution\n",
    "\n",
    "class Turing:\n",
    "    \"\"\"Represents a 1-D Turing machine.\"\"\"\n",
    "\n",
    "    def __init__(self, table, n, m=None):\n",
    "        \"\"\"Initializes the CA.\n",
    "\n",
    "        tape: map from \n",
    "        n: number of rows\n",
    "        m: number of columns\n",
    "\n",
    "        Attributes:\n",
    "        table:  rule dictionary that maps from triple to next state.\n",
    "        array:  the numpy array that contains the data.\n",
    "        next:   the index of the next empty row.\n",
    "        \"\"\"\n",
    "        self.n = n\n",
    "        self.m = n if m is None else m\n",
    "        \n",
    "        self.tape = np.zeros((n, self.m), dtype=np.int8)\n",
    "        self.head = np.zeros(n, dtype=np.int64)\n",
    "        self.head[0] = m//2\n",
    "        self.state = 'A'\n",
    "        self.table = table\n",
    "        self.next = 1\n",
    "\n",
    "    def start_string(self, init):\n",
    "        \"\"\"Starts with values from a string of 1s and 0s.\"\"\"\n",
    "        self.tape[0] = np.array([int(x) for x in init])\n",
    "\n",
    "    def loop(self, steps=1):\n",
    "        \"\"\"Executes the given number of time steps.\"\"\"\n",
    "        for i in range(steps):\n",
    "            try:\n",
    "                self.step()\n",
    "            except StopIteration:\n",
    "                break\n",
    "\n",
    "    def step(self):\n",
    "        \"\"\"Executes one time step.\"\"\"\n",
    "        if self.state == 'HALT':\n",
    "            raise StopIteration\n",
    "            \n",
    "        a = self.tape\n",
    "        i = self.next\n",
    "        a[i] = a[i-1]\n",
    "        head = self.head[i-1]\n",
    "        symbol = a[i-1, head]\n",
    "        print(symbol, self.state, end=': ')\n",
    "        new_symbol, move, self.state = self.table[symbol, self.state]\n",
    "        print(new_symbol, move, self.state)\n",
    "        \n",
    "        a[i, head] = new_symbol\n",
    "        if move == 'R':\n",
    "            head += 1\n",
    "        else:\n",
    "            head -= 1\n",
    "        self.head[i] = head\n",
    "        self.next += 1\n",
    "\n",
    "    def get_array(self, start=0, end=None):\n",
    "        \"\"\"Gets a slice of columns from the CA.\n",
    "\n",
    "        Avoids copying if possible.\n",
    "\n",
    "        start: index of first column\n",
    "        end: index of the last column plus one\n",
    "        \"\"\"\n",
    "        if start==0 and end==None:\n",
    "            return self.tape\n",
    "        else:\n",
    "            return self.tape[:, start:end]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Solution\n",
    "\n",
    "# Here's the action table for a 3-state busy beaver.\n",
    "\n",
    "table = {}\n",
    "table[0, 'A'] = 1, 'R', 'B' \n",
    "table[0, 'B'] = 1, 'L', 'A' \n",
    "table[0, 'C'] = 1, 'L', 'B'\n",
    "table[1, 'A'] = 1, 'L', 'C' \n",
    "table[1, 'B'] = 1, 'R', 'B' \n",
    "table[1, 'C'] = 1, 'R', 'HALT'"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Solution\n",
    "\n",
    "# Make the Turning machine and run it\n",
    "\n",
    "n = 20\n",
    "m = 20\n",
    "turing = Turing(table, n, m)\n",
    "\n",
    "turing.loop(n-1)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Solution\n",
    "\n",
    "# The following class draws the TM.\n",
    "\n",
    "class TuringDrawer:\n",
    "    \"\"\"Draw a Turing object using matplotlib.\"\"\"\n",
    "        \n",
    "    def draw(self, ca, start=0, end=None):\n",
    "        \"\"\"Draws the CA using pyplot.pcolor.\n",
    "        \n",
    "        start: index of the first column to be shown \n",
    "        end: index of the last column to be shown\n",
    "        \"\"\"\n",
    "        a = ca.get_array(start, end)\n",
    "        cmap = plt.get_cmap('Blues')\n",
    "        plt.imshow(a, interpolation='nearest', cmap=cmap)\n",
    "\n",
    "        # empty lists draw no ticks\n",
    "        plt.xticks([])\n",
    "        plt.yticks([])\n",
    "        \n",
    "        xs = ca.head\n",
    "        ys = np.arange(len(xs))\n",
    "        plt.plot(xs, ys, 'r.')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Solution\n",
    "\n",
    "# And here's what it looks like.\n",
    "\n",
    "td = TuringDrawer()\n",
    "td.draw(turing)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**Exercise:** This exercise asks you to implement and test several PRNGs.\n",
    "For testing, you will need to install \n",
    "`DieHarder`, which you can download from \n",
    "https://www.phy.duke.edu/~rgb/General/dieharder.php, or it\n",
    "might be available as a package for your operating system.\n",
    "\n",
    "1. Write a program that implements one of the linear congruential\n",
    "generators described at http://en.wikipedia.org/wiki/Linear_congruential_generator}.\n",
    "Test it using `DieHarder`.\n",
    "\n",
    "2. Read the documentation of Python's `random` module.\n",
    "What PRNG does it use?  Test it.\n",
    "\n",
    "3. Implement a Rule 30 CA with a few hundred cells,\n",
    "run it for as many time steps as you can in a reasonable amount\n",
    "of time, and output the center column as a sequence of bits.\n",
    "Test it.\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Solution\n",
    "\n",
    "# Here's a start, but this solution is not done.\n",
    "\n",
    "rule = 30\n",
    "m = 10001\n",
    "n = 20000\n",
    "ca = Cell1D(rule, n, m)\n",
    "ca.start_single()\n",
    "%time ca.loop(n-1)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Solution\n",
    "\n",
    "# extract the center column\n",
    "\n",
    "bits = ca.array[:, m//2]\n",
    "bits"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Solution\n",
    "\n",
    "# count the number of 0s and 1s\n",
    "\n",
    "from collections import Counter\n",
    "Counter(bits)"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.6.4"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 1
}
